


clear all
set more off
pause off
set matsize 11000



capture confirm file "Main_Results/"
if _rc mkdir "Main_Results/"


capture confirm file "Main_Results/Fig2/"
if _rc mkdir "Main_Results/Fig2/"

capture confirm file "Main_Results/Fig2/gphs/"
if _rc mkdir "Main_Results/Fig2/gphs/"


capture confirm file "Main_Results/Fig3/"
if _rc mkdir "Main_Results/Fig3/"

capture confirm file "Main_Results/Fig3/gphs/"
if _rc mkdir "Main_Results/Fig3/gphs/"



local region  Euro   Italy UK Holland Germany France  Spain 

local response R_real


local tEuro_rstar "European real interest rates (1315—2018)"
local tEuro "Response of real natural rate to a pandemic"
local tItaly "Italy"
local tUK "UK"
local tHolland "Netherlands"
local tGermany "Germany"
local tFrance "France"
local tSpain "Spain"


local tR_real "real rate"




***************************************
* LP - interest rates
***************************************

use "Data/Big_Pandemics.dta", clear
gen pandemic=1 
drop year
ren endyear year
sort year
merge 1:1 year using "Data/RealRateSince1310.dta"
sort year
replace pandemic=0 if pandemic==.
gen pandemic100k=cond(dead>=100000,pandemic,0)
cap drop Horizon
gen Horizon = _n-1 if _n<=41
tsset year



cap drop b0* se0* 
cap drop h0* l0* h1* l1*
tsset year
foreach x of local region {
	foreach y of local response {
	cap drop `x'`y'trend
	cap drop `x'`y'c
	gen b0_`x'`y'=0 if Horizon==0
	gen se0_`x'`y'=0 if Horizon==0	
	gen h0_`x'`y'=0 if Horizon==0
	gen l0_`x'`y'=0 if Horizon==0
	gen h1_`x'`y'=0 if Horizon==0
	gen l1_`x'`y'=0 if Horizon==0

	
	cap drop raw
gen raw = .
replace raw = `x'`y'
/********** State space specification starts here *************/

constraint 1	[rstar]L.rstar		= 1		/* Random Walk */	
constraint 2	[raw]rstar			= 1		/* Definitional */


cap drop zero								/* generate zero line for figures */
gen zero = 0

tsset year
cap drop r_star

/************	 State space equations	 ******************/

sspace		(rstar L.rstar, state noconstant)		/// /* r-star */
			(raw rstar, noconstant ),				/// /* raw */
			covstate(identity) constraints(1/2) difficult /*method(kdiffuse)*/ iterate(25)
			
			qui{
			/* Generate smoothed estimates of r-star */
			predict r_star if e(sample), states smethod(smooth) equation(rstar)
			
			

			cap drop Zero
			gen Zero=0
			/**** plot full sample and post-1985 for more detail ****/
			twoway   (tsline raw if year>1311, yaxis(2) lc(blue) lw(normal) ///
				 xlab(1400(200)2000) ) ///
				(tsline r_star if year>1311, yaxis(2) lc(red)  lw(thick) lp(normal) ///
				 xlab(1400(200)2000) ) ///
				 (tsline Zero if year>1311, yaxis(2) ylab(-60(20)40, labsize(large) axis(2)) lc(black)  lw(medthin) lp(normal) ///
				 xlab(1400(200)2000, labsize(large))  xtitle("Year",  size(large)) ytitle("Percent", axis(2)  size(large)) ) ///
				, title(" ") ///
				legend(col(2) order(1 2) region(color(none))  size(large)  ///
				label(1 "raw real interest rate") label(2 "real natural rate") ring(0) position(6) ) scheme(s1color) ///
				name(`x'_rstar,replace)  nodraw scale(1.2)
				
	
			cap drop `x'`y'trend
			gen	`x'`y'trend = r_star

			}


								}
								
						}
						
						
graph combine Euro_rstar, ///
	rows(1) cols(1)  ysize(2.5) xsize(4) imargin(tiny)  iscale(0.75) 
	
graph save   Main_Results/Fig2/gphs/Fig2a.gph, replace
graph export Main_Results/Fig2/Fig2a.pdf, replace

			
			
			
local region   Euro
foreach y of local response {
	foreach x of local region {
	eststo clear
		forv h=0/40 {
		
						cap drop  d`y'`h'
						gen d`y'`h' = f`h'.`x'`y'trend- l.`x'`y'trend
						reg d`y'`h' pandemic100k l(0/10).pandemic100k  l(2/11).`x'`y'trend 
												
						eststo
						replace b0_`x'`y'  = _b[pandemic100k]                      if Horizon==`h'
						replace se0_`x'`y' = _se[pandemic100k]                     if Horizon==`h'
						replace h0_`x'`y'  = _b[pandemic100k]+2*_se[pandemic100k]  if Horizon==`h'
						replace l0_`x'`y'  = _b[pandemic100k]-2*_se[pandemic100k]  if Horizon==`h'
						replace h1_`x'`y'  = _b[pandemic100k]+1*_se[pandemic100k]  if Horizon==`h'
						replace l1_`x'`y'  = _b[pandemic100k]-1*_se[pandemic100k]  if Horizon==`h'
						}
						
		esttab , keep(pandemic*)
		cap drop Zero
		gen Zero=0
		twoway ///
		(  rarea h0_`x'`y' l0_`x'`y' Horizon if Horizon~=. , fc(blue)  fi(10) lcolor(white)  lpattern(none)  )  ///
		(  rarea h1_`x'`y' l1_`x'`y' Horizon if Horizon~=. , fc(blue)  fi(20) lcolor(white) lw(none)  lpattern(none)  )  ///
		(  line  Zero    Horizon if Horizon~=. , lc(gray)  lw(med) lp(dash)                     )  ///		
		(  line  b0_`x'`y'     Horizon if Horizon~=. , lc(blue)  lw(thick) lp(solid)  ///
				xlab(0(10)40, labsize(large)) ylab(, angle(0) labsize(large)  nogrid) xtitle("Years since pandemic event ended",  size(large)) ///
				ytitle("Percent",  size(large))  ///
				legend( off ) )  ///
		, title()    ///
		xsize(2) ysize(2) scale(1.1)  scheme(s1color)  ///
		name(`x'`y',replace)  
		
graph save    Main_Results/Fig2/gphs/Fig2b.gph, replace
graph export  Main_Results/Fig2/Fig2b.pdf, replace


* Print differnet dimension of the same graph for plotting in Fig 3a
twoway ///
		(  rarea h0_`x'`y' l0_`x'`y' Horizon if Horizon~=. , fc(blue)  fi(10) lcolor(white)  lpattern(none)  )  ///
		(  rarea h1_`x'`y' l1_`x'`y' Horizon if Horizon~=. , fc(blue)  fi(20) lcolor(white) lw(none)  lpattern(none)  )  ///
		(  line  Zero    Horizon if Horizon~=. , lc(gray)  lw(med) lp(dash)                     )  ///		
		(  line  b0_`x'`y'     Horizon if Horizon~=. , lc(blue)  lw(thick) lp(solid)  ///
				xlab(0(10)40, labsize(large)) ylab(-2(1)2, angle(0) labsize(large)  nogrid) xtitle("Years since pandemic event ended", size(large)) ///
				ytitle("Percent",  size(large))   ///
				legend( off ) )  ///
		, title()    ///
		xsize(5) ysize(4) scale(1.1) ysc(r(-2.5 2.1))   scheme(s1color)  ///
		name(`x'`y'_comp,replace) 
		
graph save   Main_Results/Fig3/gphs/Fig3a.gph, replace
graph export Main_Results/Fig3/Fig3a.pdf, replace
}
}








	
	




local region   Italy UK Holland Germany France  Spain 

foreach y of local response {
	foreach x of local region {
	eststo clear
		forv h=0/40 {
		
						cap drop  d`y'`h'
						gen d`y'`h' = f`h'.`x'`y'trend- l.`x'`y'trend
						reg d`y'`h'  l(0/10).pandemic100k  l(2/11).`x'`y'trend 
												
						eststo
						replace b0_`x'`y'  = _b[pandemic100k]                      if Horizon==`h'
						replace se0_`x'`y' = _se[pandemic100k]                     if Horizon==`h'
						replace h0_`x'`y'  = _b[pandemic100k]+2*_se[pandemic100k]  if Horizon==`h'
						replace l0_`x'`y'  = _b[pandemic100k]-2*_se[pandemic100k]  if Horizon==`h'
						replace h1_`x'`y'  = _b[pandemic100k]+1*_se[pandemic100k]  if Horizon==`h'
						replace l1_`x'`y'  = _b[pandemic100k]-1*_se[pandemic100k]  if Horizon==`h'
						}
						
		esttab , keep(pandemic*)
		cap drop Zero
		gen Zero=0
		twoway ///
		(  rarea h0_`x'`y' l0_`x'`y' Horizon if Horizon~=. , fc(blue)  fi(10) lcolor(white)  lpattern(none)  )  ///
		(  rarea h1_`x'`y' l1_`x'`y' Horizon if Horizon~=. , fc(blue)  fi(20) lcolor(white) lw(none)  lpattern(none)  )  ///
		(  line  Zero    Horizon if Horizon~=. , lc(gray)  lw(med) lp(dash)                     )  ///		
		(  line  b0_`x'`y'     Horizon if Horizon~=. , lc(blue)  lw(thick) lp(solid)  ///
				xlab(0(10)40, labsize(large)) ylab(, angle(0) labsize(large)  nogrid) xtitle("Years since pandemic event ended") ///
				ytitle("Percent", size(large)) scheme(s1color) ///
				legend( off ) )  ///
		, title("`t`x''")    ///
		name(`x'`y',replace) nodraw
		

}
}





set more off
graph combine FranceR_real GermanyR_real ItalyR_real HollandR_real SpainR_real UKR_real , ///
	title()  ///
	rows(2) cols(3) imargin(tiny) ysize(7) xsize(12)  iscale(0.7) scheme(s1color) ///
	name(country_end, replace)
	
	
graph save   Main_Results/Fig2/gphs/Fig2c.gph, replace
graph export Main_Results/Fig2/Fig2c.pdf, replace
